
#include <assert.h>
#include <fstream>
#include "SAMatch.h"



//////////////////////////////////////////////////////////////////////////////
//SAMatch
void SAMatch::add(const uint3 value, const char saOffset) {
	_sa.push_back(st_SA(value.x, value.y, value.z, saOffset));
}


void SAMatch::clear() {
	_sa.clear();
}

void SAMatch::getPos(vector<unsigned int> &pos, const st_SA &sa) {
	pos.clear();
	for(unsigned int i = sa.y; i <= sa.z; i += 1) {
		pos.push_back(_saIndex[i]);
	}
}

///////////////////////////////////////////////////////////////////////////////
//CPUSAMatch

bool CPUSAMatch::isMatched(const char* seq, const char* ref, 
			const int readLength, const int numMismatch) {

	int numError = 0;
	for(int i = 0; i < readLength; i++) {
		if(seq[i] != ref[i]) {
			numError++;
		}
	}

	if(numError == numMismatch) {
		return true;
	} else {
		return false;
	}
}

int CPUSAMatch::matchReadToRef(const char* read, 
				vector<unsigned int>& pos, const char saOffset,
				const int readLength, const char* ref, 
				const int numMismatch) {
	int numMatch = 0;
	
	vector<unsigned int>::iterator it;
	for(it = pos.begin(); it < pos.end(); it++) {
		if(isMatched(read, ref + (*it) - saOffset, readLength, numMismatch)) {
			numMatch++;
		} 	
	}

	return numMatch;
}

void CPUSAMatch::match(QuerySet &querySet, Reference &ref, const int numMismatch) {
	vector<st_SA>::iterator it;
	unsigned int qid;
	vector<unsigned int> pos;
	char saOffset;
	char* read = NULL;
	const int readLength = querySet.readLength;
	int numAlignment = 0;
	cout << "numMismatch = " << numMismatch << endl;


	ofstream ofs("smallsa.txt");
	for(int i = 0; i < _sa.size(); i++) {
		ofs << _sa[i].x << " " << _sa[i].y << " " << _sa[i].z << endl;
	}
	ofs.close();

	int c = 0;
	for(it = _sa.begin(); it < _sa.end(); it++) {
		getPos(pos, *it);
		saOffset = (*it).saOffset;
		qid = (*it).x;
		read = const_cast<char*>(querySet.getRead(qid));

		numAlignment += matchReadToRef(read, pos, saOffset, 
						readLength, 
						ref.base, numMismatch);
		c++;
	}

	cout << "#sa = " << c << endl;
	cout << "CPUSAMatch::match::numAlignment = " << numAlignment << endl;
}

//////////////////////////////////////////////////////////////////////////////
//GPUSAMatch

void GPUSAMatch::match(QuerySet &querySet, Reference &ref, const int numMismatch) {
	//extract necessary data to the GPU
	const int numSA = _sa.size();
	st_SA* d_sa = NULL;
	unsigned int* d_qid = NULL;
	unsigned int* d_offset = NULL;
	unsigned int* d_length = NULL;
	GPUMALLOC((void**)&d_sa, sizeof(st_SA)*numSA);
	TOGPU(d_sa, &(_sa.front()), sizeof(st_SA)*numSA);
	GPUMALLOC((void**)&d_qid, sizeof(unsigned int)*numSA);
	GPUMALLOC((void**)&d_offset, sizeof(unsigned int)*numSA);
	GPUMALLOC((void**)&d_length, sizeof(unsigned int)*numSA);
	const int posLength = extractSA(d_qid, d_offset, d_length, d_sa);

	//expand the SA to the reference positions, modify the memory layout
	//as saIndex is too large, we have to do the SA conversion on the CPU
	unsigned int* h_pos;

	
	//cleanup
	GPUFREE(d_sa);	
}


